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Abstract 

Estimation of stochastic processes evolving in a random environment is of crucial im¬ 
portance for example to predict aircraft trajectories evolving in an unknown atmosphere. 
For fixed parameter, interacting particle systems are a convenient way to approximate such 
stochastic process. But the second level of uncertainty provided by the environment pa¬ 
rameters leads us to also consider interacting particles on the parameter space. This novel 
algorithm is described in this paper. It allows to approximate both a random environ¬ 
ment and a stochastic process evolving in this environment, given noisy observations of the 
process. It is a sequential algorithm that generalizes island particle models including a pa¬ 
rameter. It is referred by us as labeled island particle algorithm. We prove the convergence 
of the labeled island particle algorithm and we establish bound as well as time uniform 
bound for the asymptotic error introduced by this double level of approximation. Finally, 
we illustrate these results on a filtering problem where one learns a dynamical parameter 
through noisy observations of a stochastic process influenced by the parameter. 


Introduction 

This paper deals with the estimation of stochastic processes whose evolution is influenced by a 
random environment. This question is at stake in different areas. In economy, when one wants 
to estimate the option price with an unknown volatility [T] using the Black-Scholes model, one 
can consider that the option price has its evolution influenced by an unknown environment, the 
market volatility. In biology, when one wants to estimate the number of bacteria whereas the 
environment factors are unknown j2], one can model the evolution of the bacteria number as 
a stochastic process whose evolution is influenced by unknown external factors. In air traffic 
management, this modelization can also be used when one wants to predict aircraft trajectories 
evolving in an unknown atmosphere. Indeed if pilot intents and some aircraft parameters are 
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not known, actual wind and temperature evolve locally and are not perfectly known neither. 
Those atmospheric parameters which appear in the dynamic equations of the aircraft have a 
great importance to predict the future position of the aircraft. They are thus both uncertain. 
Therefore, in order to improve the trajectory prediction, one has to learn aircraft parameters 
but also atmospheric ones. It has been shown in |3] that it can be done using mode-S radar 
observations and this specific model. 

When the stochastic process evolving in the random media is Gaussian and its evolution is 
linear, the double estimation can be made using interacting Kalman filters (IKF) |4l[5]. However 
when the dynamics are non-linear, as for aircraft dynamics, an analytic resolution is not possible. 
A method based on interacting particle systems, which takes into account the randomness due 
to the environment and also the randomness coming from the process itself, was proposed by Del 
Moral in |1]. This idealized algorithm would be a sequential Monte Carlo (SMC) algorithm on 
the couple defined by the random environment and the conditional law of the process evolving in 
this random environment given the history of the environment. However, the calculation of the 
previous conditional law is not tractable in practice when the dynamics are non linear. Therefore 
another approximation level is necessary in order to estimate this conditional law. We propose 
in this paper to use interacting systems of interacting particles. These interacting systems can 
be viewed as a two-level interacting particle system. The top level particles are composed of an 
environment proposition and an empirical measure which gives an approximation of the process 
law evolving in the proposed environment. The empirical measure is obtained by the second level 
of interacting particles. This nested structure was also presented in for mean field processes. 

This algorithm can be seen as a generalization of interacting island particle models where 
each island is associated with a random parameter. Those island particle models have been 
introduced in and their statistical properties studied in [8], but without parameters. The 
first paper deals with the parallelization of interacting particle systems, the second one studies 
the asymptotic properties of the ensuing estimator. Concerning filtering problems, Chopin et al. 
in [5] introduced a kind of island particle models where each island is identified by a parameter 
proposition. They proposed an algorithm called SMC^ which is a practical version of the idealized 
iterated batch importance sampling (IBIS) algorithm introduced by Chopin in |1()| for exploring 
a sequence of parameter posterior distributions. The considered parameter did not have any 
proper dynamic whereas in the present paper the stochastic process evolution scheme depends 
on a dynamic parameter. Moreover, in the SMC^ algorithm, islands of particles grow continuously 
with time as particles ancestral lines are required to estimate the likelihood increments, and by 
their product to estimate the total likelihood. The algorithm introduced by Crisan et al. in 
m is a different version of the SMC^ which allows also the estimation of fixed parameters of a 
state-space dynamic system using sequential Monte Carlo methods. However, unlike the SMC^ 
method, the proposed algorithm by Crisan et al. operates in a purely sequential and recursive 
manner. In particular, the scheme for the rejuvenation of the particles in the parameter space is 
simpler, given that it does not need the simulation of the auxiliary particle filter from initial time 
to evaluate the likelihood. Therefore the algorithm we propose in this paper is similar to the 
algorithm of m in the sense that it is sequential in time and structured as a nested interacting 
particle filters, but different as it deals with dynamic parameters. 

In this article, we present a novel algorithm for estimating both a random environment and 
a process whose evolution depends on this environment, and study the asymptotic properties 
of the ensuing estimators. This study is of great importance to justify the convergence of this 
algorithm and also a challenging issue as it deals with error in distribution space. Therefore as 
a first step we establish bound for the asymptotic error introduced by this double level of 
approximation at every time step. As a matter of fact, the shape of the bound was predicted by 
Baehr in his thesis |5]. Then we obtain a time uniform bound for the error. From there we 
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deduce the almost sure convergence of the estimator towards the target measure. Afterwards, we 
compare the labeled island particle algorithm and interacting Kalman filters (IKF) on a filtering 
example dealing with the evolution of a mobile on a random media. In particular, it appears 
that the labeled island particle algorithm gives a better estimate of the position and the speed of 
the mobile than IKF. Finally, the labeled island particle algorithm is applied to another filtering 
problem where one learns a dynamical parameter through observations of a stochastic process 
influenced by the parameter. The theoretical results of this paper are illustrated on this example. 

Formalization of the problem through Feynman-Kac measures is given in [Section T| then 
the labeled island particle algorithm is described in [Section 2] bounds of this algorithm are 
established in [Section 3 Finally, convergence of the labeled island particle algorithm and some 
results proved in Section ^ are illustrated in [Section 5 [on two filtering examples. 


1 Feynman-Kac models in random media 

In this section, we first present an example which motivates our study, and then we introduce 
notations and models. 

1.1 Example of process evolution in random media 

In this article, one always consider stochastic processes whose evolution are influenced by their 
surrounding environment. When the environment is unknown, one can be interested in estimating 
both the environment and the law of the stochastic process itself using observations of the last 
one. Take a really simple example : a mobile evolving in whose dynamics is influenced by an 
unknown exterior force. This problem can be modeled by the following system of equations 

f Xn+l = Xn + Vn^ ^ ^ At + 0„+i At + 

\ 14+1 = Vn+BX, 

where Xn+i denotes the position of the mobile which depends on 0n+i and a Gaussian 
noise. The proper speed 14 of the mobile, is known up to a Gaussian white noise B^- a is the 
course track parameter of the mobile. The vector 0n+i is random and represents the unknown 
force acting on the position of the mobile. 

We are interested in the estimation of the state of the mobile, which depends on the parameter 
0„+i. We thus need to learn both the force, the speed and the position of the mobile. Gonsider 
now that noisy observations 14 of the mobile’s state are available. One has to estimate the 
quantity E [(Ag, Vq, 0o), • ■ •, (A„, 14, 0„) | Ig, ■ ■ • j W] • Therefore, we need to use a model which 
can tackle this issue. To this end, the formalism of Feynman-Kac models in random media is 
well adapted. In [Section 1.3[ we recall the definitions attached to this model and some important 
results. For a more detailed review see |1]. 

1.2 Notations 

Let us define some notations used in this paper. For (m, n) G 1? such that m < n we denote 
I'm, n] A {to, m-l- 1, ..., n} C Z. We will use the vector notation Um-.n — (dm, ■ ■ ■, in)- Moreover, 
IR+ and IRlj_ denote the sets of nonnegative and positive real numbers respectively, and N* the 
set of positive integers. 

N(^, E) denotes a multivariate Gaussian distribution with mean ^ and covariance matrix E. 
In the sequel we assume that all random variables are defined on a common probability space 
(n,A, P). For some given measurable space (E,£) we denote by M(E) and 'P{E) C M(E) the set 
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of measures and probability measures on (E,£), respectively. In addition, we denote by F(E) 
the set of real-valued measurable functions on (E,f) and by .^^(E) C F(E) the set of bounded 
such functions. For any v G M(E) and / e F(E) we denote by lyf = f f{x) i'{dx) the Lebesgue 
integral of / under ly whenever this is well-defined. Now, given also some other (Y, 3^) measurable 
space, an unnormalized transition kernel K from (E,£i) to (Y,3^) is a mapping from E x 3^ to 
M such that for all A € 3^, a: i—>■ K{x, A) is a nonnegative measurable function on E and for all 
a: € E, A mapstoK{x, A) is a measure on (Y,3^). If K{x,y) — 1 for all a; € E, then K is called a 
transition kernel (or simply a kernel). The kernel K induces two integral operators, one acting 
on functions and the other on measures. More specifically, let / € F(E) and n G M(E) and define 
the measurable function 

Kf-.E^x^ [ f{y)K{x,dy), 


and the measure 


vK : 3^ 9 A !—)■ / K{x, A) ^(da;). 


whenever these quantities are well-defined. Finally, let K be as above and let L be another 
unnormalized transition kernels from (Y,3^) to some third measurable space {Z,Z); then we 
define the product of K and L as the unnormalized transition kernel 

KL : E X Z B {x, A) i~G j K{x, dj/) L{y, A), 
from (E,£) to (Z,Z) whenever this is well-defined. 


1.3 Introduction of Feynman-Kac models 

Let 0„ be a random process on E® which influences the evolution of another random process 
Xn on E^. In order to avoid any confusion, all the quantities which refer to the random process 
0„ (resp. Xn) may be identified by the exponent 0 (resp. AT). Let the couple (0„,Ar„)„gN 
be a E„ = (E®, E^)- valued Markov chain of elementary transition matrix T„ form E„_i to E„ 
defined by 

Tn Xn-i),d{0n, Xn)) = M® (6l„_i, d6<„)M^ „(a;„_i, dx„), 

where and ^ are the transition kernels of the 0„ and Xn processes from E®_j^ to E® 
and from En_i to En respectively. 

Its initial distribution is given by 

rio{di9o,xo)) = rioiddo)v0„{dxo), 

with rjQ G V(Eq) and rjf^ G 'P{Eq), denoting respectively the initial distributions of 0o and Xq 
given 00 = ^'o- 

Let (G'„)„gN be a collection of bounded measurable functions from E„ to ]0, oo[. We define the 
Feynman-Kac measure associated to the couple (GmTn) with initial distribution rjQ by 


}vo,n idiiOo,xo),...,{0 n ■) ))) 

2 (' n-l 




k p=0 J 

with the normalizing constant Zn, given by 

fn—1 




[]Gp(0p,Xp) 

.p—0 


> 0 , 
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and the two path probabilities 

(d(0o,... ,0n)) = r?o®(d0o)Mf (0o,d0i)... M® d0„), 

and 

(d(a;o, ■ • ■,Xn)) = r 7 ^(da;o)Mg^i(a:o,da:i)... da:„). 

As one may have noticed, given Qo-.n = Og-.n, the sequence A„ is also a Markov chain of transition 
kernels and initial distribution 77^. Then one can associate to it another Feynman- 

Kac path measure which is called quenched. 


Definition 1.1. The quenched Feynman-Kac path measure associated to the realization Qo-.n = 
0o-.n Is defined by 




(d(a;o,.. .,Xn)) = 



(d(xo, . . ■,Xn)) ■ 


where the quenched normalizing constant 2,^ ^ is given hy 


'0:n,rL 

n—1 




Y[Gpi9p,Xp) 

.p—0 


> 0 . 


In the rest of the paper the quenched potential functions are denoted by Gg^^p and defined as 

Gg^^p • Xp G Ep I y Vsp^pixp') — Gp{0pjXpf (3) 

To get further into the dynamic, one can define the time marginal of the quenched Feynman-Kac 
measure also called the quenched Feynman-Kac distribution. 

Definition 1.2. For every realization Qo-.n = do-.m the quenched Feynman-Kac distribution flow 
'0eo:„,n 0” ^n is defined for every fn € Sh(E^) by 




X 


lifn) = 7^.„,n(/n)/7£„,„(l) 




fn{Xn)U;ZoGe,AXp) 


The distribution of A„ depends on the trajectory do-.n which is emphasized by denoting the 
unnormalized quenched Feynman-Kac distribution by 7^ „. An important result taken from 
[0, Proposition 2.6.2] is recalled below. 


Proposition 1.1. The quenched distribution sequence {rj^^ n)n6N satisfies the non linear equa¬ 
tion : 

V^0-.ri+l,n-ei ~ ^^„,ni'Tgo,„,n)^^r,+ l,n-\-li ( 4 ) 

where the mapping ^ : P(E^) — V{EAi) given by 


^ en,nAso,„,n)AXn) — x (n \Gg„AXn)'ngo,^,nAXn)- 
‘0eo.„,n\Ge,^,n) 


( 5 ) 


Defining the mapping by 




(E® X E ®+0 X P(E^) 
0n-ei)j 


'^e„,nAeo..„,n) ^K+un+l 


( 6 ) 
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The non linear recursion Q can he reformulated as 

Veo:r. + l,n+l = ^n+1 (i^n, 0n+l), ■ (V 

Remind that, for a fixed value of the random process Oo-.n, the probability measures 
can be approximated recursively thanks to an interacting particle system which 
evolves successively according to selection step with potentials Gg^^n defined in @ and transition 
kernels See [1] for further details. Now, consider that the random environment 0o:n, where 

the stochastic process evolves, is not known. Then we focus our interest on the estimation 
of the couple 

e E„ ^ (E® X V{E^)), (8) 

made up of the environment and the law of the process evolving in this environment. The tricky 
part will be to deal with the probability measure space. First, notice that, as it has been shown 
in ID, the pair process is a Markov chain. 

Proposition 1.2 (U, Proposition 2.6.3). Xn is a Markov chain with transition kernel 
defined for every function f^ G Bb(En) and {u,r]) G E„ by 

Mn(fn)iu,v)= [ M®(u,du)7„(u,$;f((u,u),77)) 

and with initial distribution fj^ G ’P(Eo) defined by 

770(d(u, - rioidu)6^x{diy). 

To this Markov chain, one may associate the Feynman-Kac distribution flow 77 ^ defined for 
every /„ e Bb(En) by 

Vuifn) =7n{7n)/7nW (9) 


where 


Inifn) =^f!o 


n—1 

fni^n) Gp{Xp) 

P=0 


and the functions Gp are non negative functions defined as follows : 

Gp : Ep —>■ [0, oo[ 

(u, 77) 1-^ Gp{u, 77) = / X Gp{u, x)T]{dx) = Lx Gu,pix)r]{dx) = f]{Gu,p)- 

P P 


( 10 ) 


Proposition 1.3 (|4], p. 86 ). For all tt G N, the sequence 77 ^ satisfies the following non linear 
recursive equation : 

7?„+i = «'„(?7 „)M„+i = $„+i(77„), (11) 

where for every fi G P(E„), the application : ’P(E„) —P(E„), is defined by 

^n(M)(7j=MGn7„)//i(G„), (V7„ e S,(E„)), (12) 


and the operator 


is defined by 

®n+i : 'P{En) 


p, I—)• 


nE„+il 

«'„(/r)M„+i. 


In the non linear case, (111 cannot be solved analytically. Therefore, in the next section, we 


introduce an interacting particle system to approximate recursively the sequence of Feynman-Kac 
probability measures ( 77 „)„gN. 
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2 Algorithm derivation 


This section is about the algorithm associated with the Feynman-Kac distribution flow 77 ^ deflned 
in (|^. One considers t he process Xn a ssociated with the pair {Gn,Mn), where the transition 


kernel M„ is deflned in 


Proposition 1.2 and the potential function G„ is deflned in (lOl. 


2.1 Idealized interacting particle model 

Let Ni be some positive integer. A Vi-interacting particle system associated with the sequence 
((G„, M„))„gN and the initial distribution rjg, is a sequence of non-homogeneous Markov chain, 

denoted by ^ , taking value in the product space E„ ^, 


A 


= = (X 




pM A F 

t •—-n - *—71 


.. X 


Ni times 


The initial state of the Markov chain consists in Ni independent random variables with 

common distribution rjQ. The interacting particle system explores the state space E„ 

and with the dynamic given to it, empirically samples the law fj^. Each particle i of the system 
consists in a random variable € E„. Therefore, the empirical process is 

deflned by 






Ni 


Ni 


(13) 


The elementary transition of the Markov chain from E„ ^ to E^+i given for any — 


i=l 

[Ni] 


G E^' by 


Ni 


G 


A^i A^i ^ / -^3 \ 

M„+i(X^,dx(j+i), thanks to 


Thus, the evolution of the particle swarm consists in two steps : a selection and a mutation. In 
the selection step, the particles (Xlj^j:^ are selected multinomially with probabil ity propo rtional 
to their potentials (G„(X^))^\. Selected particles are identified with a hat on Figure 1 Then 


the mutation step is performed independently using the kernel M^+i- The evolution scheme 
of the particles is illustrated on [Figure l] Using this algorithm one can empirically sample the 
measure 77 ^ at each time step n. Several results are available to qualify the subsequent estimator. 
However, as one may have noticed, for each 0^ the measure 77 ^ ^ corresponds to the quenched 


distributions deflned in (1.2|. That means that one should have the exact quenched measure 
associated with the parameter realization 9^.^ to use that standard particle algorithm. This can 
happen in two special cases. 

Firstly, one special case is when the transition kernel ^ is Gaussian and the initial 
distribution 77 ^ is Gaussian. Indeed, it turns out that this particle algorithm corresponds to 
the interacting Kalman Alters (IKF) (see H], |1]). That is a Xi-interacting particle model which 
is composed of Xi particles where the measure value part are Gaussian distributions. In other 
words, for each particle one iterative step of the Kalman Alter is run to update the measure. 
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Figure 1: Evolution scheme of the interacting particle system for exact measures. 


i.t. one prediction step and one correction step. Those filters are then competing through the 
selection step using the transformation defined in (12 1 . For example let us consider the case 
where ©„ is a E®-process with initial distribution 77 ® and elementary transition kernel M®. For 
a realization 0o:n of Oo:n, consider that (X„,Y„) is a IRP+'J-Markov chain, for positive integers 
{p,q), defined through the linear Gaussian system : 


f ^0,1,n ^n—1 ^0n,n “f ^ 71^1 

^ Yji — Xji “f , Tl ^ 0. 


{Ag^^n, Bg^^njCe^^n, Dg^^n) and {ag^^n, cg^^n) are respectively matrices and deterministic vectors 
of appropriate dimension which may depend on a parameter 6n- The sequences and are 
two independent white noises, independent from the initial condition Xq. There are Gaussian 
random variables whose mean and variance are given by 

Xq ~ N(TOeg_o, Sflo.o), £n ~ N(0, S^), and N(0, E^). 

In this framework, 77 ^ ^ corresponds to the conditional law of Xn given the observations 
Yo-.n-i = yo:n-i and the history of the parameter 9o-.n, also called optimal predictor. One wants 
to estimate recursively the law of the couple ( 0 „, 77 ^ ^ using observations lo:n-i = Uo-.n-i- 
For that purpose, one needs to introduce the optimal filtering which is the conditional law of 
Xn given the observations Y^-n = j/o:n and the history of the parameter 0 o:n- It turns out that 
these previous distributions are Gaussian respectively denoted by 77 ^ ^ = N{mg^^n,Yg^^n) and 
N{mg^^n,Yg^^n)- ThuS, 


‘Alg^^n — I 

Yg^^n — ^0n,n){_Xn ] 

Yg^_^^ ^n-\-l — 1^00:71.+! [(^^+1 ^n+1) (^n+1 '^g^+i ] ■ 

Moreover, the mapping ^n+i defined in which is used to update the measure valued part 
77 ^ n corresponds to a complete step of the Kalman filter evolution between predictors. This 
means that 0„_|_i), N(r 7 iti^ „)) is also a Gaussian distribution whose mean and 

covariance matrix are obtained recursively through two steps: 

N(77r,„,„,E,„,„) N(7f7,„.„, E,„.„) N(7r7,„,,,„+i, E,„,„)- 


The first one is a correction step which is given by 

f ^071,'n ~ ^0mn A Xg^^n{Yn {,Cg^^n'^07i,n A Cg^^n)) 
\ Yg^^n = {I — Kg^^nCg^^n)Yg^^n 
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where / is the identity matrix and is the classical gain matrix 

Ke„,n = ^e„,niC9„,nf {Ce„,n^e^,u{Ce„.nf + De„,n^liDg^,nfy' ■ 

The second step is the predicting step : 

J '^9n + l,n+l ^6n+l,n+l'^9n,n “t“ CL9n+l ,n+l 

\ n+1 — ^6n+i,n+l'^6j^,n{^6n+i,n+l) “t" ^n+1 ^n+1 ^'^+1) ' 

Then all the Kalman filters attached to each realization for i G |1, A^i] interact through 
their potential Gn+i(0lg.i, defined in (101 by 

= N ^ri+l ) 

where is the likelihood function defined for every a;„+i € by 

dN(Gej^ „+iX„+i, 

^n+l) 


Gg^ n+l(2;n+l) 


dN(0,E^+i) 


One finally ends up with the following expression: 


Gn+l(6'„+l, 

i,n+l)^ + ^n+l) 

= dN(0,S^+i) 


■ (14) 


See U for further details. The interacting Kalman filter for this general example is given by 
Algorithm 
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Data: rj^, (^'p)^=o. ^ 6^,0 and 

Result: Interacting Kalman approximation of 77 ^ 

/* Initialization 

for * ^ 1 to 7Vi do 

I Sample = (6'S,?7^_o) ^ ^O’ ^0 and 77 ^ ^ = N(? 7 ie, q, Eg, q) ; 

end 

for p ^ 0 to n — 1 do 

/* Selection of Kalman filters 

Sample Ip = {Ip)^i according to a multinomial distribution with probability 


Ni 


proportional to [Gp{9^,r]fk p)j^ ^ given by Q ; 

for 7 -f- 1 to Xi do 

/* Updating step for each Kalman filter 

m ji = m ji + K (¥„ — C ji m ji ) 

e/,p e/.,p s/,p^ e/,p 

E C ) E 

9/,p \ 9/,p 9/,pJ e/,p 

/* Mutation of each island 

Sample independently according to Mp_^i{9p ^,.) ; 

/* Prediction step for each Kalman filter 

TTlgi = ^4, 


>+i.p+i”^g^;_p + “e;+i,p+i 

p+i = \+„P+i^g!i P+i)^ + ^e* p+iS^+i(-Be;+,.p+i)' 

C'p ,IJ 


-p+i.p+i 
P+1 


end 

p <— p + 1 


*/ 


*/ 


*/ 


*/ 


*/ 


end 


Algorithm 1: Interacting Kalman Filter - IKF 
Secondly, when the non linear [Equation 6 ] can be solved analytically z.e. when one has access 
to the exact measure one can apply a simple interacting particle model as described in 

[Figure f| where each particle corresponds to the pair: parameter and exact measure. 

However, in most cases, this equation cannot be solved analytically, so that an additional ap¬ 
proximation is needed in order to estimate the measure ^ for each i G |1, Ai]. The next 
subsection is dedicated to the derivation of an algorithm to deal with this constraint. 


2.2 Labeled island particle model 

To tackle the case where p^ ^ analytically known, the idea consists in 

using a particle estimation of p^ ^ inside the previous interacting particle model. The ensuing 
algorithm will be called labeled island particle model in reference to the island particle model 
developed in [7], even if in the present case, each island i have a label 0 ^ whose evolution is given 
by the Markov kernel M®. The labeled island particle model consists in associating to each term 
of the sequence a sub A 2 -interacting particle system. We call sub A 2 -interacting particle 

system associated with the sequence {{Ggi „,M^ „))neN and the initial distribution p^ p, the 

sequence of non-homogeneous Markov chain taking value in the product space £ 7 ^X 2 ^ 

that is : 


ei,[N-i\ A A (ciP pz 

Sn VSn Jj — 1 VSn 5 • ■ • 5 Sn J ^ 




N 2 times 
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The initial state of the Markov chain consists in sampling N 2 independent random 

variables with common distribution 77 ^ 

The interacting particle system, denoted by explore the state space and with the 

dynamic given to it, empirically sample the law 77 ^ 

Denoting the empirical measure 


X,N2 a 
= 


N 2 






(15) 


the elementary transition of the process from to is given for any = 

{xl,...,x^^) G by 


N 2 

= n using (§ 


i=i 

N 2 N 2 


n E by (H). 

1 = 1 fc=l Af=l ) 


Define the mapping by 


: E^i X E® X 7^(E^_i) ^ ViE^) 


{{u,v),u) ^ ]\^li^n{{u,v),v){Axi), 


D =1 


then 




X,N2 




= $ 


Jf 

n+1 



(16) 


So, the evolution of the particle swarm consists in two steps: a selection and a mutation. 

In the selection step, the particles are selected multinomially with probability proportional to 
their potentials [Ggi Then the mutation step is performed independently using the 

kernel Hence, at each iteration n G N, the empirical measure approximates 

rj^i ^ when N 2 tends to 00 . Replacing 77 ^ ^ by inside the first algorithm presented, one 

gets a nested particle model named labeled island particle model. 

In order to derive precisely this algorithm, first introduce the following sequence on E„ = 
E® X ViE^), defined by ^be couple environment and empirical measure 

of the process conditionally on 0 o:n, where ~ 


Proposition 2.1. Xn is a En-Markov chain with transition kernel Mn defined for every function 
/„ G Bb{En) and {u,v) G E„ by 


MMu){u,v)= / M^{u,6.v)f^{v,^^{{u,v),v)), 


'Eg> 


(17) 


where is defined in (161, and with initial distribution fjo G P(Eo) given by 

77 o(d(u, v)) = 77 ® (du)^ x,iV 2 (di^). 


’?8o.O 
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Proof. Let a{XQ,... ,X„) stands for the cr-algebra generated by the random variables Xp, 0 < 
p<n. For all /„ e Bb{E„): 

Epo[7n{Xn)\<j(Xo,...,X^-l)] 

= [7„(0„, I a(lo,. ■., 

= E,)o[7„(0„,l>7((0n-i,0«),»7eo!!"_i,n-i) I cr{Xo,...,Xri-i)] by 
Recalling that Xn-i = n-i)> conclude that 

E^o[7„(^n) k(^0,...,^n-l)] 

' 7n(0n,$7((0n-l,0n),r?^’;!^_„„_i)M®(0„_i,d0„). 




□ 


To the Markov chain Xn, one may associate the Feynman-Kac distribution defined for every 

7 „ e Bb{En) by 

j?n(/„) = 7n(/„)/7«(l), (18) 

where 7„ is defined such that 


Ufn) = E. 


VO 


n— 1 

J^{X^)l[Gp{Xp) 

p—0 


with Gp defined in (10). 

In a similar way to fj„, the mea sure fjn satisfies a recursive equation fjn = '^n-i(jin-i)Mn, with 
the application defined in Proposition 1.3 This non linear equation can be rewritten as 


fjn = 

where the mapping <i>„ is defined as follows : 

l>„:lP(E„_i) ^ P(E„) 

r] ^ ^'„_i(77)M„. 


(19) 


( 20 ) 


As in ISection 2.1[ when this equation cannot be solved analytically one may use a particle 
model to approximate the probability measure fjn- In this case, the particles {X^ = (0)^,7;^’^^), 

i G |1, A^i|}, would be testing points on the state space E„, for (A^i, iV 2 ) G (N*)^. These particles 
explore the state space E„ and their dynamics empirically sample the law fjn when Ni gets large. 
An interacting particle system associated with the couple {Gn,Mn) and the initial distribution 
fjo, is a sequence of non-homogeneous Markov chain, x7^^\ taking value in the product space 
E^\ defined by 

= (^;)£i = ixl...,xl7^)GEn'- 

The initial state of the Markov chain consists in Ni independent random variables with 

common distribution 770 • Denote by fj^^ the empirical measure at time n, which is defined by 


~Ni A J_ X _ 


Ni 




( 21 ) 


2=1 
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The elementary probability transition, is given for any S ^n+i by 


Ni 




i=l 


The particle evolution is summarized in Figure 2 where by definitions (101 and (15), 

. N2 Afs 

GniK) = E E 

^i=i ^i=i 

The ensuing algorithm is described in Algorithm 


X, 


[JVi] _ 




Selection 

»=l 



Vi 


X,N2 \ 
1- ) 


Figure 2: Evolution scheme of the labeled island particle model. 


13 







Data: fjo, (Mp)”=o and («'p)p^o 
Result: Particle approximation of rjn 

/* Initialization */ 


for * ^ 1 to do 


Sample Xq = ido,r]fi’Q^) ~ that is 


qi „0 


so 


i.i.d X 

- \,0 


and 


N2 


— ATj Z) <5; 


t = l 


^o’‘ 


end 

for p ^ 0 to n — 1 do 

/* Selection of islands */ 

Sample Ip = {Ip)^i according to a multinomial distribution with probability 
proportional to f Z Gp (6»*, ) j ; 

for z ^ 1 to 7Vi do 

/* Selection of pcirticles inside each island */ 

Sample according to a multinomial distribution with probability 

/ i i ■ \ ^2 

proportional to ( Gp(9p’’ ,^p^’^) j ; 

/* Mutation of each island */ 

Sample independently 9p_^_-^ according to .) ; 

for j •(— 1 to N 2 do 

/* Mutation of particles */ 

Sample CSi according to p+iid"”'’" ^ ■) '■> 

end 

end 

p <— p + 1 

end 

Algorithm 2: Labeled island particle algorithm 
For every n > 0, is an estimator of pn, obtained through the labeled island particle 
model, i.e. for every /„ G Sf,(E„), 


VnHfn) 


1 _ 

^ \ ^ f (f\i X,N2 

lyr 2^ J n\^n^ 


i=l 


) 


converges to ?7n(/n) when Ni —)■ +c». 

3 bounds 

We are interested in this section in the bounds of the difference between the estimator 
and the measure rj^. To get these bounds we will use several notations. We define them before 
going further. 
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3.1 Notations 


Let (E, £) be a measurable space. For a real-valued measurable function h € we denote the 

oscillator norm osc{h) = sup(^_^/)gE 2 \h{x) — h{x')\, and Osci(E) the convex set of f-measurable 
functions with oscillations less than one. The sup norm of h is noted ||/i||oo — sup^gE \h{x)\ and 
the L^-norm ||.||p. Bi{E) C Bb{E) refers to the set of functions whose sup norm is less than 
one. For two probability measures {n,r]) € V{E)^, the Zolotarev semi-norm Ij.jly attached to a 
countable collection of bounded measurable functions in .Bi(E) is defined by 

\\^ - jj\\^ ^ snp\fj,{f) - r]{f)\. 

/eS 

To measure the size of a given class one considers the covering numbers 

defined as the minimal number of LP(/i)-balls of radius e > 0 needed to cover 
Let and /(U) denote respectively the uniform covering numbers and entropy integral 

given by 

7V'(e,5')= sup A/'(e,5',L^(/i)) (22) 

f^eviE) 

f Vlog(l+Af(s,l?))de. (23) 

Jo 

Let A denote the minimum operator and V denote the maximum operator. For a kernel M 
defined on E, the Dobrushin coefficient of M is 

/3(M) = sup osc(M(/)). 

/eOsci(E) 


Let {d{n))n>o be a sequence defined for every m > 0 by 

I d{2m)^^ = (2to)^2-™ 

1 d{2m + 1 ) 2™+1 A { 2 ni+l}^^-m+l /2 

{ y/m+l/2 

where for any positive integers {p,q) € (N*)^, {q -\-p)p = {q +p)\/q\. 

For n € N, introduce the Feynman-Kac semi-groups (resp. „) such that for all 

{xn,Xn+i) G E„ X E„+i (resp. (a;„,a;„+i) G E^ x 


Q (^n)H7^-|_i (Xyj, dx^^j-l-l), 

^resp. jbs^n-i-i) = Ge„,n(^rn)Hg^^^ dxn+i)^. 

For every {p, n) G (N)^ such that p < n, set 

Qp,n ~ Qp+l ■ ■ ■ Qni and Pp,n = Qp^n/Qp,ni^)j 

...Q^„_,^„,„andPg^^p„(/„) 4 Q^^,„,p,„(/n)/Q^^,„,p,„(l)) 

and set the normalizing constant 

Gp,n = <9p,„(l), (resp.Gep,„,p,„ = Q^^^_p_„(l)). 

Finally, set 

— A Gpn{Xp) I ^ Gq ,^^p^n{Xp)\ 

9p,n= sup_ =—1 resp. 5 ep,„,p,„ = sup 


(5?p.y„)e(Ep)2 Gp^niUp) 


Gp,i/p)e(EX)2 Gep,„,p,n{yp) I 
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In order to study the difference between and rj^, we use several results taken from |3]. 
Then, we will always assume that for all n S N, the potential functions defined in ^ 

satisfy the following condition (Gg): 

there exists a sequence of strictly positive number en{Gg) G (0,1] such that for any {xn,yn) G 

■■ 

Ge^^n{Xn) > f-n{Gg)Gg^^n{yn) > 0 (Gg) 

Therefore, for all n S N, the potential functions Gn satisfy the following condition (G): 
there exists a sequence of strictly positive number en(G) G (0,1] such that for any (xn, y„) G (E„)^ 


Gn(Xn) > ^n{G)Gn{yn) > 0 (G) 

Moreover we always assume that the collection of distributions .))_ are absolutely 

continuous with one another. That is for every n > 0 and (xn,yn) ^ (En)^, one has 


^n+l(Xnj ■) ^n+l(jjrn ■)■ 


In addition, we assume that the collection of distributions \M^ abso¬ 
lutely continuous with one another. That is for every n > 0, 9n+i G and G (E^)^, 

one has : 




Note that for time homogeneous models on finite spaces condition those conditions are met as 
soon as the Markov chain is aperiodic and irreducible. Some examples are illustrated by typical 
examples in |3]. 


3.2 bound 

Consider that for all n S N, the product space E„ = E® x 'P{E^) is equipped with the norm 
II • llg^ such that for all {u,v) G (E®)^ and (i^, ry) € (7^(E^))^, 

\\{u,r]) - {v,iy)\\^^ = \u-v\ + \\r]-iy\\:g^. 

where "Sn is a countable collection of functions in 

Theorem 3.1. For any p G N*, n G N, let /„ G Osci(E„) be a kn-Lipschitz function. As¬ 
sume that for any 9n G E®, the kernel transition ^ can he written as ^{xn-i,<ixn) = 
a;„)pe^^„(da;„) for some measurable function rnf^ ^ on EyJLi x and some prob¬ 
ability measure pg„^n G V{E^). Furthermore, assume that there exists a collection of mappings 
on E^ such that 

suPa;„_iGE^_i I logw^,„(a;„-i,a:„)| < ae„,„(a;„) 

with < oo. 

Then, the error is bounded by 


\finHfn) - Vn(.fn)\\p < + Hu)) + 2^ (24) 


where the sequence d{n) is defined in (|24[), I{'Sn) is defined in (|23|), (b(n))n>o is defined by 
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6(0) = 0 and 6(n + 1) < 5e„,nPe„+i,n+i(e '’"+1’"+^) Y. 9e,.„,q,nl3iPe,.^,q,n), 

q=0 

and a{n) is a sequence such that for all n S N*, a(n) < c [fi/2]! with c a universal constant. 
Proof. Let G Osci(E„) be a fc„-Lipschitz function, and apply triangular inequality: 

\\VuHfu)-Vn{fu)\\p < \\VnHfn) - VuH7n)\\p + \\VnHfn) - Vn(fn)\\p^ 
where ??„ is defined in (|^. Then using Theorem 7.4.4 from |5], one can bound the second term 

d{p) 


\€^Un)-VnUn)\\p < 2 ^ q,nP{P q,n) ■ 

q=0 


Therefore, in order to bound the first term, use the definitions of and in (211 and (13) 
respectively : 

\\VnHfu)-VnH7n)\\p='^rjg [\Vn" (7 n) - VnH7 




- 


Ni 

- ^ (^n, ^ ) - 7 „ (^n, } 


2=1 

N- 




p-i i/p 


p-t i/p 


As /„ is A:„-Lipschitz, it follows 

ni X,N 2 \ ~r / ni 


^ kn 

kn 


7ni0l, - fniOn, J < K {9^, - (6^, J 


X,N2 X 


where 5'„ is a countable collection of functions in Therefore one gets 




=-/v- 


P — 


— V 6- 


Ni 


i=l 


1/p 


Denoting by 9* the value at which the maximum of the Zolotarev semi-norms for i € |1, A^i] is 
reached, yields to 


ll^^'(/n) - Vn^fJllp < fcnE^o 

But, according to la, Corollary 7.4.4] (p. 247), 


1 i/p 


E. 


Vo 






where I{^n) is the entropy of defined in (23). Finally one ends up with 

\\finH7n)-V^H7n)\\p < fc„ (/(5„) + 6(n)) + 2;^ 5,,„/3(n.n)- 


□ 
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3.3 Time uniform bound 


Before stating the uniform estimate we define the following two additional conditions. 

There exists some integer to > 1 and some numbers e„(M) G ]0,1[ such that for n € N and 
_2 

(Xn,yn) G E„, one has : 

Mn,n+rniXn, ■) = Mn+1 ■ ■ ■ Mn+rn{Xn, ■) > .) 

For all i G |l,iVi], there exists some integer m-i > 1 and some numbers en{M^) G ]0,1[ such 
that for uGN, (a;„, j/„) S (E^)2, and 0l,+i,n+m, G x ... x E®+^, one has : 


M, 


X 












Theorem 3.2. Suppose that conditions ( |G|) , ({M)jn \ are met for some integer fh > 1 and some 
pair parameters {en{G),en{M)) and set e{G) = A„>oe„(G') and e{M) = An>o£„(M). 

Moreover, assume that for all i G |1, A^i] conditions {Gg I and () hold true for some se¬ 
quence of integer mi and some pair parameters (e„(Ggi), e„(M^)) and set e{Ggi) = An>oeniGgi) 
and e{M^) = An>o£niM^). Set to = ViTOi. 

Further assume that for all n > 0 and 91^ G E® the kernel transition ^ has the form 
^{xn-i,dxn) = mfi ^(Xn-i, Xn)pg^^ri(dXn) for some measurable function mfi ^ on ^ 
E^ and some probability measure pgi G P(E:^). 

Also assume that sup 2 :^_jgEX ^ jlogTO^ ^{xn-i,Xn)\ < <^6i^,n(xn) with < 00 for 

some collection of mappings agi on Eff, and set : 

< oo and pg{e^°“>) = ViPgiie^^^o'). 

n>0 


Then for any p G , any kn-Lipschitz functions G Osci(E^) one has : 


sup sup WVn^ifn) -Vnifu)\\p< 


2d(p)m 


/„60sci(E„) 

with 


k a{p) 

y]Vie(M)3e(G)2™-i 


1 + 


mpg ) 

7{Mff€{Gg)^^ 


k = snpkn, e(Mg ) = Ai£(Mgi) > 0, e(Ge) = Aie(Gei) > 0, / = sup/(S'„) < oo. 

n n>0 


Proof. 


sup_ sup _ Wvn'ifn) -Vnifn)\\p 

">0 /^gOsci(E„) 

<SUp_ SUp_ {\\fin"(fn)-Vn"(fn)\\p+\\VnHfn)-Vn(fn)\\p) 
">0 /^gOsci(E„) 


<SUp_ SUp_ WVnHfn) -Vn"ifn)\\p + Snp _ SUp _ \\v„^ (f J - Vnif n)\\p 

">0 /^gOsci(E„) ">0 /„gOsci(E„) 


From [la, Theorem 7.4.4] (p. 247), one has 
mp sup 

'‘>0 7„gOsci(E 


sup sup WVn'ifn) “ Vnifn)\\p < rTT 


2d{p)m 

^ ~ V7G£:(M)3e(G)2™-i ’ 
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since conditions ([^ and ({M)rn \ hold true. Then it follows that the only term one has to work 
on is the following. 

sup sup \\f}n^(fu)-Vn^(jr. 


/„60sci(E„) 


= sup sup 

7„eOsci(E„) 


' n)\\p 


Ni 


1 


Ni 


<—^sup_ sup_ 

i=i ^>0 /^eOsci(E„) 




As the function /„ is A;„-Lipschitz, for all i G |1, A^i], 


sup sup 

7„60sci(E„) 


= sup sup 
">0 7„gOsci(E„) 


p 

I* 


< supE^^ 

n>0 


kZ 


X,N2 


1/p 


Pi ^jP 


< sup knEjj, 


n>0 


^0 




X,N2 X 


1/p 


Set k = sup„ kn, then 


sup sup 
">0 7^gOsci(E„) 




< k sup E; 

P n>0 


Vo 


X,N2 X 


1/p 


From [la, Corollary 7.4.5] (p. 249), as one assumes that there exists > 1 for 


c F® X X F® 

^n,n-\-rai ^ '—n ^ ^ ‘-n+mj? 


sup E^ 


n>0 


Vo 


X,N2 X 


1 1/P 


<^Ei/ 


m^pg^ ) 


where I = sup/(!?„) < oo. One concludes easily. 

n>0 


□ 


4 Asymptotic analysis of the labeled island particle algo¬ 
rithm 


This section deals with the asymptotic behavior of the labeled island particle algorithm. Espe¬ 
cially, we focus on the almost sure convergence. 

Using jTheorem 3.1] obtained in jSection 3] one can easily get the al most sure con vergence of 


the double estimator f]^'^ toward under the same assumptions as in Theorem 3.1 


Theorem 4.1. Under the same assumptions as in \Theorem 3.1\ for all n > 0 and for every kn- 
Lipschitz function /„ S Osci(E„), one has 


VnHfn) ^ 7„(/n)> aS N ^OO, 


with N = A^iA ^2 such that Ni = and N 2 = ^ for all a G ]0,1[. 
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Proof. Let /„ G Osci(E„) be a fc„-Lipschitz function and e > 0 a real constant. For all p S N’ 
by Markov’s inequality, one has 


'(IVnHfJ -VnifJI > e) < 




£P 


^nH/n) -yn(fn)llp < 


Then, applying [Theorem 3.1| and noting 

C(p, n) = k„a(p)(I(J'„) + b{n)) and C{p, n) = 2d{p) 9q,nl^i.Pq,n), 

one has 

^ C{p, n) ^ C{p,n) ^ 

^ y/W j 

^ \k) 

k—0 ^ ' 

The finite sequence (sa,p(fc))fc=o defined by Sa,p{k) = {a — l/2)k + (1 — a)p/2 is bounded from 
below by 


(1 - a)p 


A “Pn 

'^a,p — 2 ^0<a<0.5 ■ ^ 


0.5<a<l 5 


so that 


VnHfn) -VnUu)\\p < 


(c{p,n) + C{p,n)^ 


'P — J^rricp 

Choose p a positive integer such that rua^p > 2 i.e. satisfying 

p > ^ if 0 < a < 0.5 
p > if 0.5 < q; < 1. 


(25) 


Hence, 


VnHfn) -VnUp 


\l< 


(C{p,n) + C{p,n)^ 

JP 


By comparison of series of non-negative general term with a convergent Riemann series, one 
concludes that the series 


ildn'ifn) - Vuifn)\ > is Convergent, 

N>0 

which implies by Borel-Cantelli’s lemma, that 

VnHJn) VnUn), asN ^OO. 


□ 

Taking such kind of N means that for a total budget N of particles, one can consider any 
decomposition (as a power of N) of the particles between islands and within each island. 
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5 Example of application 

In order to give illustration of this algorithm and of the previous theoretical results obtained, we 
present in this section two estimation problems. 

First let us recall the example of a mobile whose evolution is influenced by an unknown force 
which has been described in 0. Noisy observations of this physical systems are available. We 
resume the dynamics by the following system of equations : 


I X„+1 = “ ^At + 0„+iAt + 

] y„+i = v^ + bX 

[ Yr, = hiX„,Vn)+B^ 


(26) 


where X^+i denotes the position of the mobile in the plane, Vn the proper speed of the mobile 
and Yn their noisy observations through the observation function h, with B^ ~ N(0, S^). The 
course track of the mobile a is constant over time. The vector 0„ is a random variable and 
denotes the unknown force acting on the position of the mobile. Its equation of evolution is 
given by 


0n+i — 


0 

0 


1 

n+1 

2 

n+1 


cos 0)j 
sin0^ 


+ B 


e 

n 


with B^ ~ N(0, E®). The initial condition of the system is given by Xq N(m^g,S^Q), 
Vo ~ N(mg", ) and a = tt/2. We are interested in the estimation of the position of the mobile, 
which depends on the parameter 0„. We thus need to learn both the force, the speed and the 
position of the mobile. The tricky part is that there is no observation of the force. Here we will 
consider that the speed is a Poisson process, that is B^ is a Poisson process of intensity 0.03 
where the jumps high is given by a standard normal distribution of variance 3. Concerning B^, 
it is a Gaussian random variable such that B^ ~ N(0, E^). We present now the results obtained 
for a simulating time of 125 minutes with At = 15s. The value of the different variances are set 
to 


E® 


1 0 
0 1 


^ — ^9o,0 


1.5 0 

0 1.5 


and E^ 


0.5 0 0 \ 

0 0.5 0 

0 0 1 / 


As one can notice, to estimate the law of the couple (0„, given the observations lo:„, one 

can use Interacting Kalman filters and labeled island particle filters (LIPFs), detailed respectively 
in Algorithms [2 and 1^ We present comparative results obtained thanks to both methods. 
Concerning the labeled version, the potential of each particle is given by the density of the 
observations, that is for all Xn G and for all G E®: 


Gn{0n,Xn) OC exp 



h{Xn,Vn))’^{Y^) ^(j/„ 





On all the figures the realization of the true signal is represented by the color black, the observa¬ 
tions Y are represented by the color blue, the filtered signal obtained thanks to Algorithm with 
Ni = 100 and N 2 = 300 in red and results obtained using Algorithm in green with Ni = 100. 
On Figure 3) one realization of the signal Vn, its observed and its estimations counterparts are 
represented with respect to time. As one may observe, the true signal is well estimated by the 
technique we develop. Indeed, here the Interacting Kalman filter is not optimal as the noise 
sequence is not Gaussian. On [Figure 5| we represent the temporal evolution of the force strength 
estimation. One can notice that even if no observation is available, we are able to find back 
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the value of the true signal thanks to Algorithm whereas Algorithm retrieves only a global 
trend. Figure 4] represents the temporal evolution of one realization of the force orientation and 
its estimated counterparts. Resnlts obtained thanks to Algorithm give a better estimation 
of the trne signal than the results obtained thanks to the Algorithm From this example we 
can conclude that the labeled island particle filter is able to filter observations of the process 
while estimating the environment where the process evolves. Moreover the comparison with the 
Interacting Kalman filter algorithm shows that the labeled island particle filter is more effective 
to treat this double level estimation problem. 

Let us consider the 2-D filtering problem inspired from the growth model m- This model, which 
is a standard benchmark example in the particle filtering literatnre, is given by the following sys¬ 
tem of eqnations : 



8cos(1.2(n -I- 1)) 


-B. 


X. 




'n-l-l 


Xn + Bl 


‘n+1 


where ©o ~ N(0, cr^), Xq ~ N(0, cr^), - N(0, cr^), B^^^ - N(0, cr|) and B^ . 

We use the labeled island particle model to estimate the law of the couple (©n, Ve 
observations Yq-^, where the potential functions Gn are given by the likelihood of the observations, 
that is for all Xn S and e E®: 


N(0, CTy). 

^ given the 


Gn{0n,Xn) OC exp - 


{Yn - Xn f 

2a Y 

We present the results obtained for a simnlating time of 1000 time steps. The different variances 
are set to ag = 1, ax = 1 and ay = 10. On all the figures the realization of the true signal 
is represented in black color, the observations Y are represented in blue, and the filtered signal 
obtained thanks to Algorithmwith 7Vi = 200 and N 2 = 100 is represented in red. On Fignrej^ 
one realization of the signal © and its estimation obtained thanks to the labeled island particle 
algorithm are represented on a small period of time. As one may observe, the trne signal is well 
estimated even if no observations are available. On Fignre one realization of the process X 
is represented, its observed and its estimation connterparts. Even if the observations are really 
noisy, one is able to filter ont the noise to find back the valne of the trne signal. Indeed, as one 
may have noticed, on Figure the filtered power spectral density (in red) is closer to the black 
line, representing the “true” signal, than the observed power spectral density which has the same 
shape as a white noise for the high frequencies. Moreover, some frequencies are found even if 
there are not present in the observed signal. These two observations illustrate the convergence 
of the estimator constructed by the labeled island particle algorithm detailed in Algorithm 
Then we run 100 times the same experiment to get a sample of realizations for the true signal 
and the filtered signal. In that way one can illnstrate the theoretical resnlts obtained for the 
error bound. On fignres and 10 are presented the errors between the estimated law and 
the trne law at one time step respectively for © and X in fnnction of the nnmber of islands Ni 
and the number of particles inside each island X 2 . This error decreases both with the nnmber 
of particles and the nnmber of islands as it was snggested by the [Theorem 3.1| Concerning the 
variance of the error made between the trne law and the filtered one, on figures [m and for 
0 and X respectively, one can observe that the results obtained in [Theorem 4.1| are confirmed. 
Moreover one can notice that the variance is more influenced by the number of islands than the 
number of particles inside each island. Indeed as in [Figure 12[ the variance obtained for a fixed 
time step is varying with respect to the nnmber of islands and number of particles inside each 
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islands. But if the number of islands influences the variance, we can observe that the number of 
particles inside each island does not seem to be really influent for a given number of islands. 
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True signal 



Figure 3: Temporal evolution of the mobile’s speed, its observed and filtered counterparts. 
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Figure 4: Temporal evolution of the force orien¬ 
tation in rad 


Figure 5: Temporal evolution of the force 
strength in m.s~^ 
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Figure 6: Temporal zoom on one realization 
of the 0 process and its estimated counter¬ 
part 


Figure 7: Temporal zoom on one realization 
of the X process, its observed and estimated 
counterparts 
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Power spectral densities 
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Figure 8: True, filtered and observed power spectral densities for one realization of the process 
X over 1000 time steps. 



Figure 9: Evolution of the estimation’s er¬ 
ror for the law of 0 for 100 realizations of 
the process in function of and X 2 



Figure 10: Evolution of the estimation’s er¬ 
ror for the law of X for 100 realizations of 
the process in function of and X 2 



Figure 11: Evolution of the variance esti¬ 
mation’s error for the law of 0 for 100 re¬ 
alizations of the process in function of 
and X 2 


Figure 12: Evolution of the variance esti¬ 
mation’s error for the law of X for 100 re¬ 
alizations of the process in function of 
and X 2 
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